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We investigate particle spectra and elliptic flow coefficients in relativistic heavy ion collisions by 
taking into account the distortion of phase space distributions due to bulk viscosity at freezeout. 
We first calculate the distortion of phase space distributions in a multi-component system within 
the Grad's fourteen moment method. We find some subtle issues when one matches macroscopic 
variables with microscopic momentum distributions in a mult i- component system and develop a 
consistent procedure to uniquely determine the corrections in the phase space distributions. Next, 
we calculate particle spectra by using the Cooper- Frye formula to see the effect of the bulk viscosity. 
In spite of the relative smallness of the bulk viscosity, we find that it is likely to have a visible effect in 
particle spectra and elliptic flow coefficients. This indicates the importance of bulk viscosity together 
with shear viscosity if one wants to constrain the transport coefficients with better accuracy from 
experimental data. 

PACS numbers: 25.75.-q, 25.75.Nq, 12.38.Mh, 12.38.Qk 



I. INTRODUCTION 

One of the major discoveries at the Relativistic Heavy 
Ion Collider (RHIC) in Brookhavcn National Laboratory 
(BNL) [1] is that the quark-gluon plasma (QGP) behaves 
like a "perfect liquid" [2], which attracts much theoret- 
ical interest in various fields. Hydrodynamic modeling 
of heavy ion collisions plays an important role in de- 
ducing the above fact [3, 4]. Anisotropy of radial flow, 
namely elliptic flow [5], is found to be large in Au+Au 
and Cu+Cu collisions at RHIC energies [6-8]. Elliptic 
flow is quantified by the second harmonics of azimuthal 
angle distribution, V2 = (cos2(0 p )) [9]. One reproduces 
the centrality, pseudorapidity, and transverse momen- 
tum dependences of «2 data reasonably well by employ- 
ing the Glauber type initial conditions and implementing 
hadronic dissipative effects in ideal hydrodynamic models 
[10]. However, when one replaces initial conditions from 
the Glauber model to the ones from a model based on 
the color glass condensate, elliptic flow coefficients over- 
shoot the experimental data due to eccentricity larger 
than that from the conventional Glauber model calcula- 
tions [10] . Given a fact that the initial conditions in rel- 
ativistic heavy ion collisions are not known precisely, the 
discrepancy strongly suggests effects of viscosities in the 
QGP phase are required to reproduce the V2 data in this 
particular case. Even though ideal hydrodynamic models 
have been successful in reproduction of the «2 data, the 
next important tasks are to constrain the transport coef- 
ficients by comparison of theoretical predictions with the 
experimental data and to comprehensively understand 
the transport properties of the QGP. 

Shear viscosity has already been taken into account in 
several hydrodynamic simulations to investigate its effect 
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on elliptic flow in relativistic heavy ion collisions [11-16]. 
On the other hand, few studies are done so far to investi- 
gate the effect of bulk viscosity in dynamical simulations 
of relativistic viscous fluids [17-21]. Many years ago, it 
was found that the bulk viscous coefficient £ has a promi- 
nent peak in the vicinity of the quantum chromodynam- 
ics (QCD) phase transition due to reduction of the sound 
velocity [22] . Recently, the importance of bulk viscosity 
is realized again in the context of violation of scale invari- 
ance in QCD [23]. The smallness of the sound velocity is 
intimately related with violation of scale invariance. So 
it would be a universal feature that the bulk viscous co- 
efficient becomes large around the critical region. From 
a phcnomenological point of view, large bulk viscosity 
could trigger a violation of applicability of hydrodynamic 
framework around (pseudo-)critical temperature [24] and 
support the success of hybrid approaches [10, 25, 26] in 
which the hydrodynamic description of the QGP is fol- 
lowed by the kinetic description of the hadron gas. 

In this paper, we focus on the effects of bulk viscosity 
in the late stage of relativistic heavy ion collisions. To 
demonstrate them, we calculate transverse momentum 
spectra and differential elliptic flow coefficient. In ideal 
hydrodynamic calculations, the Cooper-Frye formula [27] 
is conventionally employed to calculate the particle spec- 
tra at freezeout. In the hybrid calculations, it plays a 
role in transition of description from macroscopic hydro- 
dynamics to microscopic hadronic cascade models. In 
both cases, the Cooper-Frye formula is indispensable in 
comparison of hydrodynamic results with experimental 
data. In dissipative hydrodynamics, viscous effects are 
taken into account in the Cooper-Frye formula in two 
ways: One is a variation of flow profile as a matter of 
viscous correction to the dynamical evolution, the other 
is distortion of phase space distribution from its equi- 
librium form. So far, no full three dimensional viscous 
hydrodynamic simulations are available. Given this situ- 
ation, we take flow profile from a full (3+l)-dimensional 
ideal hydrodynamic simulation [10, 28] assuming varia- 
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tion of it due to viscosity is not significant. We estimate 
the correction from bulk viscosity as well as from shear 
viscosity to the distribution function by matching macro- 
scopic quantities with the ones calculated in the kinetic 
theory and see how this affects transverse momentum 
(pt) spectra and transverse momentum dependence of 
elliptic flow coefficient V2(j>t)- 

This paper is organized as follows. In Sec. II, we give a 
brief overview on the rclativistic kinetic theory by putting 
emphasizes on subtle issues in its application to relativis- 
tic multi-component gases and to the ones with zero (net- 
baryon) number density limit. In Sec. Ill, we focus on a 
hadron resonance gas model, which is widely exploited in 
hydrodynamic analyses of relativistic heavy ion collisions, 
and introduce a simple model for shear and bulk viscous 
coefficients. We demonstrate the effects of bulk viscosity 
at freezeout is possibly large by calculating the prcfac- 
tors in correction of distribution functions. In Sec. IV, 
we calculate £>T-spectra and V2{pt) by taking flow pro- 
files from three cases: (1) Bjorken flow, (2) a blast wave 
model, and (3) a full (3+l)-dimensional ideal hydrody- 
namic simulation. Section V is devoted to conclusion. 



II. THE RELATIVISTIC KINETIC THEORY 

We would like to express viscous corrections to 
the phase space distribution for a relativistic multi- 
component gas in terms of macroscopic variables by 
matching the macroscopic variables with the ones ob- 
tained from the kinetic theory. This is a generalization 
of the treatment developed by Israel and Stewart [29] for 
a single component system to a multi-component system. 

We will see that some non-trivial facts arise when a 
multi-component gas is considered. Firstly, one usu- 
ally considers a scalar term, a vector term, and a trace- 
less tensor term for a momentum expansion of the non- 
equilibrium part of the distribution and assumes that the 
trace part of the tensor term can be absorbed in the scalar 
term. However, this is no longer the case for a multi- 
component system, because the equivalence of the scalar 
term and the trace part holds only in a single component 
gas. Secondly, the number of equations which determines 
the modification of the distribution does not change upon 
taking a zero net baryon density limit. This means we 
do not have to suffer the loss of conditions to obtain the 
unique expression of the distortion, nor do we need to 
resort to the quadratic ansatz [16]. Consideration of the 
limit is of practical importance since the limit is often 
employed in hydrodynamic models of the hot QCD mat- 
ter created in relativistic heavy ion collisions. 



A. Distortion of the Distribution Function 

Tensor decompositions of the energy-momentum ten- 
sor T^" and the net baryon number current Ng give the 
definitions of the macroscopic variables. The decomposed 



terms are put into the equilibrium part and the non- 
equilibrium part: 

TfiV = T v> + ST n» i (!) 

where 

Tg v = e u»u" - PA"" , (2) 
§T nv = + 2W^u^ +Tr fiU , (3) 

and 

K = N£o + SN£, (4) 

where 

K = n B0 u^ 8N£ = V". (5) 

The time-like and space-like projection operators are, re- 
spectively, defined as and A^ v = g^ v — u^u v with the 
Minkowski metric g^ v = diag(+, — , — , — ). We also define 
the dot product as a ■ b = a^b v g^ . The round brackets 
stand for symmetrization, A^B^ = \{M L B V + A V B^). 
e = UpTtf v u„, P = -^A^T^", and n B0 = u M iV£ de- 
note the energy density, the hydrostatic pressure, and 
the charge density, respectively. = A fJ, a ST al3 uis 

is the energy current, — A^SNg the charge cur- 
rent, II = -^A^ST^ the bulk pressure, and = 
r<M"> = [I(A^A^ + A^A^)-±A^A Q/3 ]T Q ' 3 the shear 
stress tensor. Generally eo and Ubo would have dissipa- 
tive counterparts Se and Sub, but such terms are not 
considered here because we employ the Landau match- 
ing conditions which demand Se = u^ST^ u v = and 
Sub = Ufj_5Ng = 0. These conditions are necessary to 
make the system thcrmodynamically stable in the first 
order theory. We will see the details later. For reviews 
on relativistic viscous hydrodynamics, see Ref. [4]. 

Kinetic theory demands for a system with multi- 
species particles 

N v _ v- f b i9i d s p 

i ' 

The label i denotes a particle species, gi and bi are the 
degeneracy and the baryon number. We define Sfi = 
fi — /P to be the deviation of the distribution function 
from its equilibrium form, /P = [exp (Ei^iEs.) - e]" 1 = 
(cxpyo — e) _1 . Here e = +1 for bosons and e = — 1 
for fermions. T and /ib denote the temperature and the 
baryon chemical potential. 

Now we expand y defined by fi = (expy — e) _1 up 
to the second order in momentum and estimate viscous 
corrections to the distribution function in the Grad's 14- 
moment method 

Sf* = -m + ef^Me^+p^e^), (8) 
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where and £ M „ are coefficients of the expansion. The 
numbers of unknown variables are 4 and 10, respectively, 
because the latter should be symmetric since any anti- 
symmetric term cancels out upon taking contraction with 

s^v is often considered to be traceless and the scalar 
term e is introduced instead [29]. The apparent lack of 
the scalar term in the expression (8) is due to the fact 
that we set the tensor term to have a non-zero trace. 
The numbers of unknown variables are the same in both 
cases. If we consider a single component system such as 
a pion gas, the trace part of the tensor term is equivalent 
to a scalar term e since the trace part can be separated 
as — -^-g^y and the metric produces a scalar p fi p u g fJ ,u — 

™2 



ties 



= E 



9id 3 p 



(permutations)] J mn , 



(10) 



= E 



(2tt) 3 £ 
V [(A MlM2 ...A M2 ™ _1 

n 

(permutations)] J mn 



(11) 



However, this is not the case for a multi-component 
system because all the viscous correction tensors e, 
and e^u are macroscopic quantities in the sense that they 
do not depend on any particular particle species. On the 
other hand, the trace part has mass dependence, which 
means it is particle-species dependent: 



8ft 



Tr(e^) 



Tr(e M1/ ) 



-/o(l + ^o) 



UfU 



Tr(£^)mf 



9V 



Pi Pi ^V-v 



(9) 



Here is a traceless tensor. The equivalence of the 
scalar term and the trace term holds for a single com- 
ponent system because we can set a variable Tr^,,) to 
be Since the same prescription does not work for 

a multi-component gas where we have the additional in- 
dex i, it is problematic whether to consider a non-zero 
trace tensor term, or the combination of a scalar term 
and a traceless tensor term in a system of multi-specics 
particles. We will see in the next section that only the 
former is relevant in the case of a 16-component hadron 
resonance gas. 

Our aim here is to uniquely determine the explicit 
forms of and £ M „. As we mentioned earlier, we in- 
troduce the Landau matching conditions u^ST^u^ = 
and Uf+SNg = 0, which ensure the thermodynamic stabil- 
ity in the first order theory (see Appendix A for further 
details). They demand the energy density and the num- 
ber density to match the equilibrium densities. Together 
with the kinetic definitions of the macroscopic variables, 
we now have 14 equations in total to determine 14 un- 
knowns in Sf. 

We first define J^iMa- ••/*«. and Jau^-m™ to express the 
conditions explicitly in terms of thermodynamic quanti- 



J mn is a coefficient of the expansion of J^to—Hm D y 
(in — 2?i) u^'s and n A Ml "s. J mn is defined in the same 
way, but has the baryon number bi as an weight factor. 
These quantities should be distinguished in a mixture 
system of baryons, anti-baryons, and mesons because 
they contribute differently to the energy-momentum ten- 
sor and the baryon number current as seen in Eqs. (6) 
and (7). 

The Landau matching condition for the energy- 
momentum tensor is simplified by using the expressions 
defined above: 



^E 



9id 3 p 



PtPiSfi 



(2^) 3 ^ 

-J 30 £* - (ho - Ju)£** - J4iTr(e). 



(12) 



From now on we employ the notations e* = £** = 

e„ v u^u u , Tr(e) = e£, and A^e ut = A» v u a e ua . The 
other conditions can be expressed in a similar fashion: 





n 



(Jso - J 3 i)£** - JsiMe), (13) 



■hi — 2^42 )e* 



— J20£ 

- J 3 iWe v - 2 JuWe„, 
- J 2 iA» v s u - 2 j 31 A» v e u *, 
-2J i2 e { ^. 



-J 42 Tr( £ ), (14) 

(15) 
(16) 
(17) 



The conditions (12)-(17) are classified into three indepen- 
dent sets of equations: (a) the definition of the bulk pres- 
sure and the Landau matching conditions for the scalars 
£*, £**, and Tr(e), (b) the definitions of the energy cur- 
rent and the charge current for the vectors A pi/ £„ and 
A^e,,*, and (c) the definition of the shear stress ten- 
sor for the tensor e^ u K The equations are immediately 
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solved for each group: 



B. Discussion 



D U, 



Tr(e) = e£ = B 3 TL, 



(nv) 



(18) 
(19) 
(20) 
(21) 
(22) 
(23) 



where the prefactors Di and Bi in Eqs. (18)-(23) are 
functions of J mn and J mn - Thus e M and e M „ are deter- 
mined to be 



= D n.u» + D 1 w li + b x v IJt , 

+ A (ttI/ (Tr(e)-£**)/3 



(24) 



+ 2u (m A 1/)q £ q * + £(,») 



= (-B A MI , + Boit^ii^)!! 

+ 2Bi« (M TV 1/) +2B lU{tl V v) +fl 2 7r^, (25) 

where J5o = (B3 — Bq)/3 is employed. The prefactors are 
explicitly expressed as 



#0 


= 3(J40 J31 - 




(26) 


A) 


= (J30^30 — 


JioJ2o)J 3 , 


(27) 


A) 


= 3(J4i J20 - 


_ ^30^3l)J3 j 


(28) 




= — 2J3iJ7 2 


1 

; 


(29) 




= 2J±iJ 2 , 




(30) 




= J21J2 i 




(31) 


Si 


= —J31J2 1 


) 


(32) 


S 2 


= ^f 1 , 




(33) 



where 

J3, = 5 J30 J42J30 + 3 J31 J40 J31 + 3 J41 J41 J20 

— 3J31 J41 J30 — 3J30 J4i>/3i — 5J40J42J20, (34) 

.72 = 2J31 J31 — 2J41 J21, (35) 

Ji = -2J 42 . (36) 

The prefactor for the shear stress tensor (33) has the 
same form as shown in Ref. [29], whereas those for the 
bulk pressure (26)-(28) are different from the ones in the 
reference because the finite trace tensor term is consid- 
ered instead of the scalar term. Of course, those reduce 
to the same formulas in a single component gas. 

It is worthwhile to mention that, if we define the heat 
current g' 1 as 

q , = — = w„- ^±±V M , (37) 

J 21 n B 

Eqs. (21) and (22) can be expressed as Diq u — V^/ubT 
and Biq^, respectively. 



Some comments are in order here. 

Firstly, we give a consideration to the special case 
where the Landau frame, i.e., W 1 = and the zero net 
baryon density limit Ng — ► are employed. These are 
often taken for hydrodynamic analyses of the QGP in 
relativistic heavy ion collisions because the numbers of 
baryons and anti-baryons are roughly the same for such 
events. In our framework, they correspond to estimating 
only the contributions of the shear viscosity and the bulk 
viscosity. 

Apparently the Landau condition of the charge cur- 
rent (13) vanishes upon taking the limit because J rnn 
does. This could induce an ambiguity because the num- 
ber of conditions for the determination of the prefactors 
decreases. But if we expand J mn around the baryon 
chemical potential fis = 0, we have 



Jmn = + Hb 



dJ„ 



dfi B 

Then Eq. (13) is reduced to 
<9(J 3 o - J31) 



0(p B ). 



(38) 



fj, B =0 



dJ 2 , 



din 



cV 3 i 
dj-iB 



Tr(e) 



0, (39) 



which remains finite in the zero baryon density limit. 
This enables us to determine the prefactors for the bulk 
pressure uniquely even after the limit is taken. 

Secondly, we discuss the validity of the quadratic 
ansatz in the distortion of the distribution. If one makes 
the quadratic ansatz, i.e., £ M „ = C\-k^ v + C^A^II [16], 
then it violates the Landau conditions because 



- J 4 iTr( e ) + 0, 

- J 3 iTr(e) ± 0. 



(40) 
(41) 



Therefore, one should not make the ansatz if one em- 
ploys the Landau matching conditions in hydrodynamic 
equations. Note that, in the quadratic ansatz, £»*'s in 
Eqs. (12) and (13) vanish since is perpendicular to 
it' 1 in this ansatz. Also the correction from bulk viscos- 
ity is naively expected to be smaller than that of shear 
viscosity, because C 2 = —\C\ holds. 

Thirdly, it is important to explicitly consider a multi- 
component system for estimations of the viscous correc- 
tions Eqs. (26)-(33) because the deviation of the distri- 
bution for the i-th particle species Sfi generally depends 
on whether it is the only component or one of the com- 
ponents in a multi-component gas. This comes from the 
fact that J mrl 's, and thus the prefactors, include informa- 
tion of all particle species like thermodynamic variables. 
It should be noted that the difference becomes negligible 
if one treats the first order theory and assume that the 
shear viscosity 77 is proportional to the entropy density s, 
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because in the Boltzmann approximation we have 



<*/shear = B 2 X 2? ? V <A1 U l/) X rfrf 
1 



2J42 
1 

2J42 
= C\T). 



2,1 



42 



T 3 



(42) 



Here we used the relation J42 ~ sT 3 in the Boltzmann 
approximation [29]. We see that C"(T) is a function de- 
pendent on the index % and the temperature but inde- 
pendent of what other components are in the gas. 



III. THE MODEL 




0.02 0.04 0.06 0.08 0.1 0.120.14 0.160.18 0.2 

T (GeV) 



A. Equation of State and Transport Coefficients 

Following the discussion in the previous section, we 
estimate the thermodynamic quantities and the prcfac- 
tors appearing in Sf. As the model of the equation of 
state, we consider the 16-component hadron resonance 
gas, which has mesons and baryons with mass up to 
A(1232). Hereafter we take the Landau frame and set 
the system to be baryon free. The models for the shear 
viscosity r\ and the bulk viscosity £ are taken from Refs. 
[30] and [31]: 



1 

c 



1 






s, 


4- 















3 



(43) 
(44) 



where 



' jr— is the sound velocity. The factor a 



in the bulk viscosity is set to 15 unless there is further 
notification. Figure 1 shows temperature dependence of 
rj/s and £/s for a hadronic resonance gas model. In the 
temperature range relevant to relativistic heavy ion col- 
lisions 0.1 < T< 0.2 GeV, (/s is several factors smaller 
than rj/s. So one might expect that the effect of the bulk 
viscosity is small. This is not necessarily true since the 
corrections in Sf appear as a combination of transport 
coefficients and the pref actors discussed later. 

Figure 2 shows the numerical results of the prcfactors 
-Do, -Bo, Bo, and B2 appearing in Sf as functions of tem- 
perature T in a hadronic resonance gas model. Back 
to the discussion whether we should consider a non-zero 
trace tensor term or a trace part of a scalar term, we find 
that only the former is relevant for the 16-component 
hadron resonance gas case. This is because, if we choose 
the combination of the scalar term and the traceless ten- 
sor term, -Do, -Bo, and -Bo diverge due to a change of the 
sign in the denominator (34) at temperature below 0.2 
GeV. 

The numerical code to generate the prefactors in e M 
and £ M „ will become available [32] . One can set the com- 



FIG. 1: (Color online) rj/s (solid) and £/s with a = 15 (dot- 
ted) as functions of temperature in a resonance gas model. 



ponents of a resonance gas and calculate the prcfactors 
as well as any moment of /o(l + e/o), J-mn- 



B. Particle Spectra and Flow 

In hydrodynamic analyses of the QGP created in rel- 
ativistic heavy ion collisions, the Cooper-Frye formula 
is employed at frcezeout. This converts the macroscopic 
variables into microscopic distributions, which enables us 
to compare the hydrodynamic results with experimental 
data. From a phcnomcnological point of view, it works 
as an interface from a hydrodynamic model to a cascade 
model. £>T-spectra of the i-th particle is given as, using 
the Cooper-Frye formula [27], 



d 2 N r 



d 2 p T dy (27r) 3 J s 



(45) 



where px and y denote the transverse momentum and 
the rapidity. The equilibrium distribution /o is given as 



fo 



1 



,(pi-U-bill B )/T 



(46) 



The elliptic flow parameter is the coefficient of the sec- 
ond harmonics cos(2<^ p ) in a Fourier expansion of the 
azimuthal momentum distribution, where 4> p is the az- 
imuthal angle in momentum space. It is defined as 

( x S d yJS*WpCoa(2<l> p ) J* 

J ay J Q "-(Pp d<t> p p T dp T dy 

Viscous corrections are taken into account via (a) vari- 
ation of the flow and (b) modification of the distribu- 
tion function. Because we do not have a full (3+1)- 
dimensional viscous hydrodynamic flow, we focus on the 
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FIG. 2: (Color online) (a) The prefactor in e M for the bulk pressure, Do, as a function of temperature, (b) The prefactors in 
e M „ for the bulk pressure, Bo (dashed), and Bo (dash-dotted), and that for the shear stress tensor B2 (solid) as functions of 
temperature. 



latter in this study. Profiles of the flow u M and the freeze- 
out hypersurface dcr^, which are necessary for calcula- 
tions of the Cooper- Frye formula (45), are taken from 
(a) the Bjorken model [33], (b) a blast wave model [34], 
and (c) a (3+l)-dimcnsional ideal hydrodynamic simula- 
tion [28]. 

IV. RESULTS 

We estimate the particle spectra and the elliptic flow 
coefficient of the negative pion with the mass m = 0.139 
GeV. The freezeout temperature Tf is set to 0.160 GeV. 
This temperature is sufficiently near the QCD (pseudo- 
critical temperature. At Tf = 0.160 GeV, 77 = 1.31 x 
10~ 3 GeV 3 and C = 4.37 x 10~ 4 GeV 3 when a = 15. In 
this study the Navier-Stokes limit is taken for the shear 
stress tensor and the bulk pressure, which means 7r Miy = 
2r ? V (At iW and II = -(d^ where V^ 1 = A^d u . 



It is noteworthy that the Cooper-Frye formula for this 
model can be analytically expressed in the case of the 
Boltzmann approximation. See Appendix B for the de- 
tails. 

Figure 3 shows the particle spectra with corrections 
from the bulk viscosity or the shear viscosity. The mean 
transverse momentum, (pt), is lowered by the bulk vis- 
cosity in the case of the Bjorken flow. This can be inter- 
preted from the sign of the shear stress tensor and the 
bulk pressure as follows. The bulk pressure works as a 
negative pressure, i.e., the pressure in the energy mo- 
mentum tensor is effectively reduced to P — \U\ because 
the system is expanding in the longitudinal direction and 
expansion scalar <9 M u M = 1/r is positive. As a result, the 
number of particles with lower momenta increases. In 
a similar way, the shear viscosity is naively expected to 
enhance (pr) in the mid-rapidity region because the pres- 
sure in the radial direction is increased by the decrease 
of the pressure in the longitudinal direction due to the 
fact that the shear stress tensor is traceless. 



A. Bjorken Flow 

We first consider the Bjorken model. We employ the 
expanding coordinates (r,r,(f>,ri s ), where proper time r, 
radius 7', azimuthal angle <f>, and space-time rapidity rj s 
are defined in the relations t = t cosh 77,,, x = rcostfi, 
y = 7'sin0, and z = Tsmhrj s . In this frame, the Bjorken 
flow is written as 

u T = 1, u r = 1/ = vP" = 0. (48) 

The radius of the nuclei and the freezeout time are set to 
i?o = 10.0 fm and r = 7.5 fm, respectively Elements of 
freezeout hypersurface only have the radial component: 

da T = rdr] s rdrd(f>, da r = da^ = da Vs = 0. (49) 



B. Flow from Blast Wave Model 

As a blast wave model, we follow Ref. [34] 

u r = u ^-[l + W2Cos(2</))]e( J Ro-r) (50) 

U T = y/1 + (U r ) 2 (51) 

= u n ° = (52) 

where uq = 0.55 and ui = 0.2. The radius of the nuclei 
and the freezeout time are set to Ro = 7.5 fm and t = 
5.25 fm. The profile of the freezeout hypersurface is the 
same as Eq. (49) in the Bjorken model. 



7 




~i , i i i i i i , i , i i i , i i , i i i i i i i i , i , i . i . . . i . \ 
0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 
P T (GeV) 

FIG. 3: (Color online) Viscous corrections to pr-spectra with 
Bjorken flow. Solid, thick dashed, and thick dotted lines are, 
respectively, the results without any viscous corrections, with 
the effect of bulk viscosity (a = 15), and with the effect of 
the shear viscosity. Thin dashed and dotted lines show that 
the absolute value of the ratio of the correction to the ideal 
spectrum becomes greater than 0.5. 

Again the bulk viscosity lowers the spectrum as shown 
in Fig. 4. The elliptic flow parameter v-zipr), on the 
other hand, is enhanced by the bulk viscosity. The 
counter-intuitive result can be interpreted as a fact 
that the slope of differential i>2(pt) can be given as 
dv2{pr) I dpT — V2/ {pt) [35]. The average V2 is not much 
affected by the viscous correction in the case of the 
blast wave model. Effects of the shear viscosity can be 
explained likewise. 



C. (3+l)-Dimensional Ideal Hydrodynamic Flow 

In the last example, we demonstrate the effect of 
the bulk viscosity by employing flow profiles from a 
full (3+l)-dimensional ideal hydrodynamic simulation in 
Au+Au collisions at y/SNN = 200 GeV. Impact parame- 
ter is taken to be 7.2 fm. For details of this particular hy- 
drodynamic model, see Refs. [10, 28]. Here the hadronic 
equation of state in this hydrodynamic model is exactly 
the same as the one considered here, i.e., a hadronic res- 
onance gas model up to the mass of A(1232). Note that 
these flow profiles together with temperature distribu- 
tions are publicly available [36]. It should be also noted 
that we have checked that H/P and ir^ v /T^ u at each 
freezeout position have moderate values, at most of the 
order of unity. 

We see in Fig. 5 that both shear and bulk effects lower 
the spectrum in the high px region. The enhancement 
of V2{pt) by the bulk viscosity is due to the decrease in 
mean px of the spectrum as discussed in the previous 



subsection. The non-triviality lies in the correction of 
the shear viscosity: It lowers the pT-spectrum and still 
suppresses i^Cpt)- This behavior was also observed for a 
(2+l)-dimensional viscous hydrodynamic flow [12]. 

The viscous corrections in these calculations may have 
been overestimated for two reasons. Firstly, we neglected 
Su^ = — Uq by considering the ideal hydrodynamic 
flow. Viscosity tends to make the system moderate 
and consequently thermodynamic forces, in particular, 
derivatives of the flow are expected to become smaller 
in the case of viscous hydrodynamic flow. This would 
reduce the amount of corrections considered here. Sec- 
ondly, the bulk pressure and the shear stress tensor are 
estimated in the first order theory, which naively means 
that they correspond to asymptotic values at the time 
much longer than relaxation times in the second order 
theory. 

Figure 6 shows comparison of p-r-spectra and V2(pr) 
between with and without the quadratic ansatz. Neither 
the amount of the correction nor its pt dependence seems 
to be similar to each other. The quadratic ansatz greatly 
underestimates the effects of the bulk viscosity. There- 
fore, the proper procedure to obtain the prefactors in 5f 
should be made as discussed in Sec. II to correctly cal- 
culate the distortion of particle spectra due to the bulk 
viscosity. 

Finally when both shear viscosity and bulk viscosity 
arc considered, the slope of the particle spectra becomes 
steeper but that of V2{pt) becomes moderate compared 
with that of the ideal distribution. The reason for the 
latter would be that the effect of the bulk viscosity and 
that of the shear viscosity accidentally cancels each other 
in the low px region. 



V. CONCLUSIONS 

We estimated the viscous corrections on the phase 
space distribution of a rclativistic gas in a multi- 
component system. We found that generalization to a 
multi-component gas involves some subtleties. Firstly, 
the trace part of the tensor term and the scalar term 
in the distortion of the distribution function are equiv- 
alent for a single component gas, but not for a multi- 
component one. Numerical calculations suggest that 
one should have the trace part in the case of the 16- 
componcnt hadron resonance gas. Secondly, if one take 
zero net baryon density limit, one also has to take the 
limit for the Landau matching condition for the baryon 
number current because it does yield a finite relation. 
Since the number of equations does not change, we can 
uniquely determine the prefactors for the bulk pressure 
appearing in the distortion of the distribution. It is not 
desirable to introduce an additional ansatz that violates 
the matching conditions, since such an ansatz makes the 
system thermodynamically unstable in the first order the- 
ory. Thirdly, the deviation of the distribution for the i- 
th gas component dfi is generally different depending on 
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P T (GeV) P T (GeV) 

FIG. 4: (Color online) Viscous corrections to (a) pr-spectra and (b) V2{pr) in the blast wave model. Solid, thick dashed, and 
thick dotted lines are, respectively, the results without any viscous corrections, with the effect of bulk viscosity (a = 15), and 
with the effect of shear viscosity. Thin dashed line shows that the absolute value of the ratio of the correction to the ideal 
spectrum becomes larger than 0.5. 




FIG. 5: (Color online) Viscous corrections to (a) pr-spectra and (b) «2(pt) with ideal hydrodynamic flow. Solid, thick dash- 
dotted, thick dashed, and thick dotted lines are, respectively, the results without any viscous corrections, with the effect of 
bulk viscosity with a = 15/2 and with a = 15, and with the effect of shear viscosity. Thin dash-dotted and dashed lines show 
that the absolute value of the ratio of the correction to the ideal spectrum becomes larger than 0.5. 



whether it is the only component or one of the compo- 
nents. This comes from the fact that the prefactors in e M 
and e^, like thermodynamic variables, include informa- 
tion of all the components in the gas. 

For hydrodynamic models of the QGP created in rel- 
ativistic heavy ion collisions, the Cooper-Frye formula is 
necessary to convert macroscopic variables into micro- 
scopic distribution at freczcout. This enables us to com- 
pare the hydrodynamic results with experimental data, 
or to see further development of the hadronic matter in 
a cascade model. Non-equilibrium effects are taken into 
account in the formula via the variation of the flow and 



the modification of the distribution. We focused only 
on the latter and estimated the viscous corrections on 
the pr-spectra and the elliptic flow «2(pt) following the 
aforementioned discussion. Profiles of the flow and the 
freezeout hypersurface are taken from the Bjorkcn model 
with cylindrical geometry, the blast wave model with az- 
imuthally anisotropic flow, and a (3+l)-dimcnsional ideal 
hydrodynamic model for the numerical calculations. We 
found that corrections of the bulk viscosity due to the 
distortion of the distribution have a visible effect on par- 
ticle spectra and elliptic flow coefficient V2(pt)- This 
implies the importance of the bulk viscosity in constrain- 
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FIG. 6: (Color online) Viscous corrections to (a) pr-spectra and (b) i>2(pt) with (thick dotted) or without (thick dashed) the 
quadratic ansatz. Solid line is the result from an ideal hydrodynamic model. Thin dashed line shows that the absolute value 
of the ratio of the correction to the ideal spectrum becomes larger than 0.5. 
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thick dash-dotted, and thick dashed lines are, respectively, the results without any viscous corrections, with the effect of both 
shear and bulk viscosities with a — 15/2, and with a = 15. Thin dash-dotted and dashed lines show that the absolute value of 
the ratio of the correction to the ideal spectrum becomes larger than 0.5. 



ing the transport coefficients with better accuracy from 
experimental data. 

Quantitatively speaking, the viscous effects might have 
been overestimated because we approximated the shear 
stress tensor and the bulk pressure as the ones in the 
first order theory and unlike in the second order theory, 
no relaxation effects are considered. As for the case of 
the estimation with an ideal hydrodynamic flow, the vis- 
cous corrections might be further exaggerated because 
the derivatives of the ideal hydrodynamic flow are gener- 
ally larger than those of real viscous hydrodynamic flow. 
Considering the fact that the shear viscosity also has non- 
trivial effects on particle spectra depending on the flow 



profile, a full (3+l)-dimcnsional viscous hydrodynamic 
flow is necessary to see more realistic behavior of px~ 
spectra and v<z(jpt)- 
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APPENDIX A: THE LANDAU MATCHING 
CONDITIONS 

The Landau matching conditions ensure the thermo- 
dynamic stability in the first order theory. If the entropy 
current s M had a term proportional to Tlu^ in the non- 
equilibrium case, the derivative 9(u A1 s AI )/9II|n=o would 
be finite. This means that the system is not in a maxi- 
mum entropy state, i.e., is not thermodynamically stable. 
The matching conditions take out such unwanted terms. 
One can explicitly show by inserting phase space distri- 
bution into the definition of the entropy current in the 
relativistic kinetic theory 



E 



E 



J^^[(l + ef) log (1 + ef) + f log f 

d3p p^ f /^- b ^ + o[(sm 



{2ir) 3 E 



T 

— «*E J J^fE? A/< 



T 



E J 5 f< 



(Al) 



Here we neglect higher order terms in the third line for 
simplicity. By using the Landau matching conditions, 
one can eliminate correction terms proportional to u M 



x" + y u " 5t 1 

„„„ 1 



up 



W» - y. B V» 

= s ir H — 

T 



(A2) 



This expression includes only non-equilibrium corrections 
at the first order with respect to dissipativc currents 
which arc perpendicular to vf*. This is exactly the gener- 
alization of the entropy current in ideal hydrodynamics 
to the one in viscous hydrodynamics as discussed in Ref. 
[29]: 



T 



(A3) 



formula for the Bjorken model in the Boltzmann approx- 
imation following Ref. [34]. The flow and the freezeout 
hypersurface of the Bjorken scaling solution with cylin- 
drical geometry are 



1, u r = u* = u 71 ' 



0. 



and 



(Bl) 



da T = Tdrj s rdrd4>, da r = da^ = da Vs = 0. (B2) 



Since the momentum p T is given by mx cosh (y — r] s ) 
where ttlt is the transverse mass, the particle spectra 
is written as 



Ro 

rdr 



d 2 N g 
d 2 p T dy (2tt) 3 J Jq 
x m T cosh (y - r] s )f 



Tdij s 



(B3) 



Here we dropped the index i of particle species for sim- 
plicity. The viscous corrections is now analytically ex- 
pressed as 



d 2 SN hulk 
d 2 p T dy 
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(2tt)3 
C 



nRl ■ Tm T 2Ki(x) 



a Do K 2 (x) + 1 
B m + -rrrnr 

2 Ki(x) 

B Q -B Q a f K 3 (x) „ 
K x {x) 



, (B4) 



and 



d 2 SN, 



hear 



(Pprdy 



(2tt) : 

2ri 
x — B 2 
3r 



-irRl ■ Tm T 2K 1 {x) 



-1 - 



Pt 



(B5) 



where x = and the modified Bcsscl function is defined 

as 



APPENDIX B: ANALYTIC EXPRESSIONS OF 
SPECTRA 



K n (x) = e- xcosh(t) cosh (nt)dt. (B6) 
Jo 



We can analytically express the corrections of both The bulk pressure and the shear stress tensor are esti- 
shear viscosity and bulk viscosity on the Coppcr-Fryc mated in the first order theory 
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